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This paper presents a brief review of the current status of ultrahigh energy cosmic ray obser- 
vations and discusses nearby starburst-like galaxies as their possible origin. 



1 Energy Spectrum, Mass Composition, and Distribution of Arrival Directions 

In the almost structureless, fast falling with energy, inclusive cosmic ray (CR) spectrum, three 
kinematic features have drawn considerable attention for a long time. These features, known 
as the knee, the ankle, and the ultraviolet cutoff, are the only ones in which the spectral index 
shows a sharper variation as a function of energy, probably signaling some "new physics". A 
plethora of astrophysical^ and exotic^ explanations have been proposed to address the pro- 
duction mechanism at the high end of the spectrum. In the absence of a single model which is 
consistent with all data, the origin of these particles remains a mystery. Clues to solve it are not 
immediately forthcoming from the data, particularly since various experiments report mutually 
inconsistent results. 

In recent years, a somewhat confused picture of the energv spectrum has been emerging. 
Since 1998, the AGASA Collaboration has consistently reported^ a continuation of the spectrum 
beyond the expected Greisen-Zatsepin-Kuzmin (GZK) cutoff^ which should arise at about 
10 1 GeV if CR sources are at cosmological distances. In contrast, the most recent results from 
HiRes^ describe a spectrum which is consistent with the expected GZK feature. This situation 
exposes the challenge posed by systematic errors (predominantly arising from uncertainties in 
hadronic interaction models^) in these types of measurements. 

Recent HiRes data have been interpreted as a change in CR composition, from heavy nuclei 
to protons, at ~ 10 GeV0 This is an order of magnitude lower in energy than the previous 
crossover deduced from the Fly's Eye dataP The end-point of the galactic flux is expected 



to be dominated by iron, as the large charge Ze of heavy nuclei reduces their Larmor radius 
(containment scales linearly with Z) and facilitates their acceleration to highest energy (again 
scaling linearly with Z) . The dominance of nuclei in the high energy region of the Galactic flux 
carries the implication that any changeover to protons represents the onset of dominance by an 
extra-galactic component. The inference from this new HiRes data is therefore that the extra- 
galactic flux is beginning to dominate the Galactic flux already at ~ 10 9 GeV. Significantly, 
this is well below Eqzk ~ 10 10 7 GeV, the threshold energy for resonant P7cmb — > A + — > Nir 
energy- loss on the cosmic microwave background (CMB), and so samples sources even at large 
redshift. 

The dominance of extra-galactic protons at lower energy can be consistent with recently 
corroborated structures in the CR spectrum. A second kne e, recognized originally in AGASA 
data^is now confirmed by the HiRes-MIA Collaboration^^ At 10 8 ' 6 GeV, the ene rgy spectrum 
steepens from E~ 3 to E~ 3,3 . This steepening at the second knee can be explained'^ by energy 
losses of extra-galactic protons over cosmic distances, due to e + e~ pair-production on the CMB. 
The theoretical threshold of the energy-loss feature occurs at 10 8 ' 6 GeV, and therefore allows 
for proton dominance even below this energy. However, the HiRes 

coincident with the energy of the second knee (from about 50% protons just below to 80% 
protons just above), and therefore argues for the beginning of extra-galactic proton dominance 
at the second knee. Another feature in the CR spectrum is the ankle at ~ 10 10 GeV where 
the spectrum flattens from E 1-3 - 3 to E~ 27 . This has been commonly identified with the onset 
of the extra-galactic flux in the past. In the aftermath of the new HiRes data, the ankle can 
now be interpreted as the minimum in the e + e~ energy-loss feature. These changes in viewing 
the onset of the extra-galactic component have spurred a refitting of the CR data down to 
10 8 ' 6 GeV with appropriate propagation functions and extra-galactic injection spectrd?^ The 
major result is that the injection spectrum is significantly steeper than the standard cx E~ 2 
predicted by Fermi engines. This appears to contradict the latest HESS measurements of the 
spectral slope of both, identified and unidentified gamma-ray sources (albeit at lower energies, 
see M. Lemoine-Gourmard in this proceedings). 

At the highest energies, the arrival directions of CRs are expected to begin to reveal their 

origins. If the CR intensity were isotropic, then one should expect a time-independent flux from 

each direction in local detector coordinates, i.e., declination and hour angle. In that case, a 

shower detected with local coordinates could have arrived with equal probability at any other 

time of a shower detection. For any point of the celestial sphere, the expected shower density 

can be estimated if the exposure in each direction can be obtained. This implies that celestial 

anisotropies can be easily discerned by comparing the observed and expected event frequencies 

at each region. Although there seems to be a remarkable agreement among experiments on pre- 
nyl 

dictions about isotropy on large scale structure,- - this is certainly not the case when considering 
the two-point correlation function on small angular scale. The AGASA Collaboration reports 
observations of eve nt cl usters which have a chance probability smaller than 1% to arise from a 
random distribution^^ whereas the recent analysis reported by the HiRes Collaboration showed 
that their data are consistent with no clustering among the highest energy eventsP" In our opin- 
ion, it is very important to rigorously define the corresponding budget of statistical significance 
and search criteria prior to studying the data, since defining them a posteriori may inadvertently 
introduce an undetermined number of "trials" and thus make it impossible to assign the correct 
statistical significance to the search result. In this direction, with the aim of avoiding accidental 
bias on the number of trials performed in selecting the angular bin, the original claim of the 
AGASA Collaboration was re-examined considering only those events observed after the original 
claim.16 This study showed that the evidence for clustering in the AGASA data set is weaker 
than was previously supposed, and is consistent with the hypothesis of isotropically distributed 
arrival directions. 



On a separate track, evidence has been presented for neutral particles (with energy ~ 
10 9 GeV) from the Cygnus spiral arm. Specifically, AGASA has revealed a correlation of the 
arrival direction of the CRs to the Galactic Plane (GP) at the 4<r level^^ The GP excess, wh ich 
is roughly 4% of the diffuse flux, is mostly concentrated in the direction of the Cygnus regional 
Evidence at the 3.2a level for GP enhancement in a similar energy range has also been reported 
by the Fly's Eye Collaboration^^ Interestingly, the full Fly's Eye data include a directional 
signal from the Cygnus region which was somewhat lost in an unsuccessful attempt to relate it 
to 7-ray emission f rom Cygnus X-3P^ The complete isotropy up to about 10 7,7 GeV revealed 
by KASCADE data^ vitiate direction -preserving photons as primaries. Therefore, the excess 
from the GP is very suggestive of neutrons as candidate primaries, because the directional sig- 
nal requires relatively-stable neutral primaries, and time-dilated neutrons can reach the Earth 
from typical Galactic distances when the neutron energy exceeds 10 9 GeV. If the Galactic mes- 
sengers are neutrons, then those with energies below 10 9 GeV will decay in flight, providing a 
flux of co smi c antineutrinos above 1 TeV that should be observable at kilometer-scale neutrino 
telescopes?^] A measurement of the P-flux will supply a strong confirmation of the GP neutron 
hypothesis. 

Another point of agreement among the experiments is the absence of a photon component at 

no 

the highest energies 1 ^ This was most recently supported by an analysis of AGASA data: abov e 
10 112 GeV, less than 50% (65%) of the primary CRs can be photons at 90%CL (95% CL)(P 

The confusing experimental situation regarding the GZK feature, mass composition, and the 
small-scale clustering in the distribution of arrival directions should be resolved by the Pierre 
Auger Observatory (PAO), which will provide not only a data set of unprecedented size, but 
also the machinery for controlling some of the more problematic systematic uncertainties. The 
first PAO site is now op erat ional in Malargiie, Argentina, and is in the process of growing to 
its final size of 3000 km 2 P^ At the time of writing, 12 telescopes and about 700 wat er tanks 
were operational. The first analyses of data from the PAO are currently underway^ Figure ^ 
shows the arrival directions of all events recorded from January to July 2004. The pixels have a 
size of 1.8 degrees and the map was smoothed with a Gaussian beam of 5 degrees. On 21 May 
2004, one of the larger events recorded by the surface array triggered 34 stations. A preliminary 
estimate yields an energy ~ 10 11 GeV and a zenith angle of about 60°. First physics results will 
be made public in the Summer of 2005 at the 29th International Cosmic Ray Conference. 



2 Possible Astrophysical Origin 



Supernova remnants (SNRs) are thought to be the main source of both CR ions and electrons 
with energies below the knee (see Reft^for a review, and references therein for details). CRs of 
low energies are also expected to be accelerated in OB associations, through turbulent motions 
and collective effects of star winds. When the Larmor radius approaches the accelerator 
size, it becomes very difficult to magnetically confine the CR to the acceleration region, and 
thus to continue the accelerating process up to higher energies. If one includes the effect of 
the characteristic velocity /3c of the magnetic scattering centers, the above argument leads to 
the general condition, E max ~ 2/3 cZeB r L , for the maximum energy acquired by a particle 
travelling in a medium with magnetic field B. In the case of one-shot acceleration scenarios, the 
maximum reachable energy turns out to have a quite similar expression to the shock acceleration 
case. In what follows, we focus on only one possible astrophysical origin of ultrahigh energy 
CRs, and refer the reader to^D for thorough discussion on other candidates. 
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Figure 1: All event (from January to July 2004) skymap in Galactic coordinates. Un its x are related to the 
number of events per pixel n according to: n = 330 x + 1.71f)P^ 

2.1 The two nearest starbursts 

Starbursts are galaxies undergoing a large-scale star formation episode. They feature strong 
infrared emission originating in the high levels of interstellar extinction, strong Hll-region-type 
emission-line spectrum (due to a large number of O and B-type stars), and considerable radio 
emission produced by recent SNRs. Typically, starburst regions are located close to the galactic 
center, in the central kiloparsec. This region alone can be orders of magnitude brighter than 
the center of normal spiral galaxies. From such an active region, a galactic-scale superwind is 
driven by the collective effect of supernovae and particular massive star winds. The enhanced 
supernova explosion rate creates a cavity of hot gas (~ 10 8 K) whose cooling time is much 
greater than the expansion time scale. Since the wind is sufficiently powerful, it can blow out 
the interstellar medium of the galaxy, preventing it from remaining trapped as a hot bubble. As 
the cavity expands, a strong shock front is formed on the contact surface with the cool interstellar 
medium. The shock velocity can reach several thousands of kilometers per second and ions like 
iron nuclei can be efficiently accelerated in this scenario, up to ultrahigh energies, by Fermi's 
mechanism 1 ^ If the super-GZK particles are heavy nuclei from outside our Galaxy, then the 
nearby (~ 3 MpJ22l) starburst galaxies M82 (I = 141°, b = 41°) and NGC 253 (I = 89°, b = -88°) 
are prime candidates for their origin. 

M82 is probably the best studied starburst galaxy, l oca ted at 3.2 Mpc. The total star 
formation rate in the central parts is at least ~ 10 M Q vr -1 ^ ! The far infrared luminosity of the 
inner region within 300 pc of the nucleus is ~ 4 x 10 10 L^P^j There are ~ 1 x 10 7 M© of ionized 
gas and ~ 2 x 10 8 M Q of neut ral gas in the IR source^^^l The total dynamical mass in this 
region is ~ (1-2) x 10 9 M P^The main observational features of the starburst can be modelled 
with a Salpeter I MF e xtending from 0.1 to 100 Mq. The age of the starburst is estimated in 
~ (1 - 3) x 10 7 yrPH Around ~ 2.5 x 10 s M (i.e. ~ 36 % of the dynamical mass) is in the 
form of new stars in the burstP^I The central region, then, can be packed with large numbers of 
early-type stars. More than 60 individual compact radio sources have been detected within the 
central 200 pc of NGC 253,^ half of which are SNRs. The supernova rate is estimated to be as 
high as 0.2 — 0.3 yr _1 , comparable to the massive star formation rate, ~ 0.1M© yr . Assuming 
that the star formation rate has been continuous in the central region for the last 10 9 yrs, and 
a Salpeter IMF for 0.08-100 M , the bolometric luminosity of NGC 253 is consistent with 1.5 



xlO Mq of young stars'^" Based on this evidence, it appears likely that there are at least tens 
of millions of young stars in the cen t r al reg ion of the starburst. These stars can contribut e to 
the 7-ray luminosity at high energies ^ '■'^ l A TeV detection was reported by CANGAROO^^, 
but has been yet unconfirmed by other experiments. A multiwavelength model of NGC 253, 
predicting less gamma-ray flux than that detected by CANGAROO (thus, predicting that HESS 
will not detect the starburst) has been presented in RefP 



2.2 Acceleration processes in starbursts 

Due to the nature of the central region, and the presence of the superwind, the escape of the iron 

nuclei from the central region of the galaxy is dominated by convection." Nuclei can then escape 

through the disk in opposite directions along the symmetry axis of the system, being the total 

path travelled substantially shorter than the mean free path. Once the nuclei escape from the 

central region of the galaxy they are injected into the galactic-scale wind and experience furt her 

acceleration at its terminal shock. CR acceleration at superwind shocks was first proposed in 40 

in the context of our own Galaxy. The scale length of this second shock is of the order of several 

tens of kpcp^ so it can be considered as locally planar for calculations. The shock velocity 

can be estimated from the empirically determined superwind kinetic energy flux E sw and the 

mass flux M generated by the starburst through: E BW = l/2Mv^ h . The shock radius can be 

approximated by r ~ v s ^t, where r is the starburst age. Since the age is about a few tens of 

million years, the maximum energy attainable in this configuration is constrained by the limited 

acceleration time arising from the finite shock's lifetime. The photon field energy density drops 

to values of the order of the CMB, and consequently, iron nuclei are safe from photodissociation 

while energy increases to ~ 10 11 GeV. 

To estimate the maximum energy that can be reached by the nuclei, consider the superwind 

terminal shock propagating in a homogeneous medium with an average magnetic field B. If we 

work in the frame where the shock is at rest, the upstream flow velocity will be vi (|vi | = v s h) and 

the downstream velocity, V2. The magnetic field turbulence is assumed to lead to isotropization 

and consequent diffusion of energetic particles which then propagate according to the standard 

transport theory^^ The acceleration time scale is then: t acc = ^ where n is the upstream 

v i 

diffusion coefficient which can be written in terms of perpendicular and parallel components to 
the magnetic field, and the angle 9 between th e (u pstream) magnetic field and the direction of 
the shock propagation, k = k\\ cos 2 9 + k± sin 2 oW** Since strong turbulence is expected from the 
shock we can take the Bohm limit for the upstream diffusion coefficient parallel to the field, i.e. 
7jE/ ZeB\, where B\ is the strength of the pre-shock magnetic field and E is the energy of 



the Z-ion. For the k± component we shall assume, following Biermannp^that the mean free path 
perpendicular to the magnetic field is independent of the energy and has the scale of the thickness 
of the shocked layer (r/3). Then, k± = 1/3 r{v\ — V2) or, in the strong shock limit, k± = rvf/12. 
The upstream time scale is i acc ~ r/(3vi), r/3v\ = 4/u 2 {E j(3ZeB\) cos 2 9 + ru 2 /12 sin 2 9) . 
Thus, using r = v\t and transforming to the observer's frame one obtains £ max ~ jZeBv^ h T ~ 

\ZeB^-T. 

2 M 

The predicted kinetic energy and mass fluxes of the starburst of NGC 253 derived from the 
measured IR luminosity are 2 x 10 42 erg s _1 and 1.2 M© yr -1 , respectively?^ The age is in the 



"The relative importance of convection and diffusion in the escape of the CRs from a region of disk scale 
height h is given by the dimensionless param eter, q = Voh/Ko, where Vo is the convection velocity and kq is 
the CR diffusion coefficient inside the starbursm^ When q < 1, the CR outflow is difussion dominated, whereas 
when q > 1 it is convection dominated. For the central region of NGC 253 a convection velocity of the order of 
the expanding SNR shells ~ 10000 km s _1 , a scale height h ~ 35 pc, and a reasonable value for the diffusion 
coefficient ko ~ 5 x 10 26 cm 2 s _1 , lead to q ~ 216. Thus, convection dominates the escape of the particles. The 
residence time of the iron nuclei in the starburst results £res ~ h/Vo ~ 1 x 10 11 s. 



range 5 x 10 7 to 1.6 x 10 8 yr (also valid for M82)*r^ F ina lly, the radio and 7-ray emission from 
NGC 253 are well matched by models with B ~ 50/x although it might be at least a factor 
of 4 bigger in the innermost region of the nucleusP^ With these figures, assuming a conservative 
age r = 50 Myr, one obtains a maximum energy for iron nuclei of E^ x > 3.4 x 10 20 eV. 



2.3 Testing the starburst hypothesis 

For an extragalactic, smooth, magnetic field of ~ 15 — 20 nG, diffusive propa gation of nuclei 
below 10 20 eV evolves to nearly complete isotropy in the CR arrival directionsP^EHl Thus, we 
could use the rates at which starbursts inject mass, metals and energy into superwinds to 
get an estimate of the CR-injection spectra. Using equal power per decade over the interval 
10 18 - 5 eV < E < 10 20 ' 6 eV - we obtain a source CR-luminosity w 3.5 e 10 53 eV/s where e 

is the efficiency of ultrahigh energy CR production by the superwind kinetic energy flux. With 
this in mind, the energy- weighted, approximately isotropic nucleus flux at 10 19 eV is given 
E3j ( E ) = {A,/dD{E) E lidt J * ~ 2 - 3 x 10 26 e/,eV 2 m- 2 s- 1 sr- 1 , where h = I M 82 + W 253- 
To estimate the diffusion coefficient we used B n Q = 15, ^Mpc = 0.5, and an average Z = 20. 
We fix e/* = 0.013, after comparing the above equations to the observed CR-flux. Note that 
the contribution of 1m82 and Ingc 253 to I* critically depends on the age of the starburst. The 
relation "starburst- age/superwind-efficiency" derived from Eq. (|2.3[) . leads to e « 10%, if both 
M82 and NGC 253 were act ive for 115 Myr. The power requirements may be reduced assuming 
contributions from M82 "B"P1 

In the non-diffusive regime (i.e., 

10 20.3 e y < E < 1Q 20.5 e y^ the accumu i ate d deflection 

angle from the direction of the source in the extragalactic -B-field is roughly 10° < 9 < 20°@^ 
The nuclei suffer additional deflection in the Galactic magnetic field. In particular, if the Galactic 
field is of the ASS type, the arrival direction of the 4 highest energy CRs can be traced backwards 
to one of the starbursts^ Figure El shows the extent to which the observed arrival directions 
of the highest energy CRs deviate from their incoming directions at the Galactic halo because 
of bending in the magnetic field. It is seen that trajectories for CR nuclei with Z > 10 can be 
further traced back to one of the starbursts, within the uncertainty of the extragalactic deviation. 

We now attempt to assess to what extent these correlations are consistent with chance 
coincidence. The deflections in the extragalactic and Galactic fields (regular and random com- 
ponents) may be assumed to add in quadrature, so that the angular sizes of the two sources are 
initially taken as cones with opening half-angles between 40° and 60°, which for the purpose of 
our numerical estimate we approximate to 50°. The global structure of the field will substan- 
tially diminishing the effective solid angle. The combined deflections in the I and b coordinates 
concentrate the effective angular size of the source to a considerably smaller solid angle. As 
a conservative estimate, we retain 25% of this cone as the effective source size. By randomly 
generating 4 CR positions in the portion of the sky accessible to the existing experiments (dec- 
lination range 5 > —10°), an expected number of random coincidences can be obtained. The 
term "coincidence" is herein used to label a synthetic CR whose position in the sky lies within 
an effective solid angle Sl c fT of either starburst. f2 e ff is characterized by a cone with opening 
half-angle reduced from 50° to 24° to account for the 75% reduction in effective source size 
due to the magnetic biasing discussed above. For the 4 observed events, with zero background, 
the Poisson signal mean 99% confidence interval is 0.82 — 12.23. Thus our observed mean for 
random events, 0.81 ± 0.01, falls at the lower edge of this interval, yielding a 1% probability for 
a chance occurrence. Assuming an extrapolation of AG AS A flux (E 3 J Q b s (E) ) up to 10 20 ' 5 eV, 
the event rate at PAO, with an aperture A ~ 7000 km 2 sr for showers with incident zenith 
angle less than 60°, is given by ^ = A E 3 J(E) % « | (E 3 J(E)) - ^] « 5.3 yr" 1 , 

where E\ = 10 20 3 eV and E2 = 10 20,5 eV. Considering a 5-year sample of 25 events and that 
for this energy range the aperture of PAO is mostly receptive to CRs from NGC253, we al- 
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Figure 2: Left: Directions in Galactic coordinates of the four highest energy CRs at the boundary of the Galactic 
halo. The diamonds represent the observed incoming directions. The circles and arrows show the directions of 
neon and iron nuclei, respectively, before deflection by the Galactic magnetic field. The solid line is the locus of 
incoming directions at the halo for other species with intermediate atomic number. The stars denote the positions 
of M82 and NGC253. The dashed lines are projections in the (I, b) coordinates of angular directions within 20° 
and 30° of the starbursts. Right: Curves of constant probabilities in the two-dimensional parameter space defined 
by the size of the cone and the minimum number of events originating within the resulting effective solid angle. 

low for different possibilities of the effective reduction of the cone size because of the Galactic 
magnetic field biasing previously discussed. In Fig. |2]we plot contours of constant probabilities 
(P = 10 -4 , 10 -5 ) in the two-dimensional parameter space of the size of the cone (as a fraction of 
the full 50° circle) and the minimum number of events originating within the resulting effective 
solid angle. The model predicts that after 5 years of operation, even if 7 are observed, it would 
rule out a random fluctuation at the 10 -5 level. 

Other galaxies, presenting a more extreme starburst behavior, like those luminous and ultra- 
luminous infrared galax ies which are also expected to be detected by space and ground-based 
gamma-ray experiments^ 1 may additionally contribute to the ultrahigh energy CR flux. 
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